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Abstract. We consider a semi-classically scaled Schrodinger equation with 
WKB initial data. We prove that in the classical limit the corresponding 
Bohmian trajectories converge (locally in measure) to the classical trajectories 
before the appearance of the first caustic. In a second step we show that 
after caustic onset this convergence in general no longer holds. In addition, we 
provide numerical simulations of the Bohmian trajectories in the semiclassical 
regime which illustrate the above results. 



1. Introduction 

1.1. WKB asymptotics. We consider the time-evolution of a quantum mechan- 
ical particle described by a wave function ip B (t, •) G L 2 (R d ; C) and governed by the 
Schrodinger equation (in dimcnsionless form): 

(1.1) K^^-^Af + VWf, V £ | t =o = V>c), 

where x £ Mr, t e M + , and V(x) € K a given potential (satisfying some regularity 
assumptions to be specified below) . In addition, we have rescaled the equation such 
that only one semi-classical parameter < e <C 1 remains. 

The classical limit of quantum mechanics is concerned with the asymptotic be- 
havior of solutions to ( 1.1 1 as e — > 0+. A possible way to describe these asymptotics 
is based on the time-dependent WKB method, where one makes the ansatz (see, 
e.g., PHHH] for more details) 

(1.2) 1> e (t,x) = a e {t,x)e iS ^' s 

for some e-independent (real-valued) phase function S(t, x) G M and a (in general 
complex valued) amplitude a e (t, x) € C satisfying 

a £ ~ a + eai + e 2 ci2 + . . . , 

in the sense of asymptotic expansions. Assuming for the moment that a e and S are 
sufficiently smooth, one can plug ( 1.2 1 into ( 1.1 1 and compare equal powers of e in 
the resulting expression. This yields a Hamilton- J acobi equation for the phase 

(1.3) d t S+^\WS\ 2 + V(x) = 0, ^=0 = 50, 
and a transport equation for the leading order amplitude 

(1.4) d t a + Va-VS + |A5 = 0, a\ t=0 = a . 
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Note that the latter can be rewritten in the form of a conservation law for the 
leading order particle density p := \a\ 2 , i.e., 

(1.5) d t p + div(pVS') = 0. 



(1.6) 



The main problem of the WKB approach is that (1.3) in general does not ad- 
mit unique smooth solutions for all times. This can be seen, from the method of 
characteristics, where one needs to solve the following Hamiltonian system 

■ X(t,y) =P(t,y), X(0,y)=y, 

P(t, y) = -VV(X(t, y)), P(0, y) = VS (y). 

By the Cauchy-Lipschitz theorem, this system of ordinary differential equations can 
be solved at least locally in-time, which yields a flow map Xt : y h-> X(t, y). If wc 
denote the corresponding inverse mapping by Y t : x i-> Y(t,x), i.e., Y t o X t = id, 
then the phase function S satisfying (1.3) is found to be (see, e.g., [TT|) 

'1, 



(1.7) S(t,x) = S (Y(t,x)) + J ^-|P(r,y)| 2 - V(X(r,y))j dr\ y=Y{t x) . 

Given such a smooth phase function S, one can, in a second step, integrate the 
amplitude equation (1.4) along the flow X t to obtain the amplitude in the form 

(1-8) a(t,x) = -= 

^J t (Y{t,x)) 

where Jt{y) '■= detV y X(t, y) is the Jacobian determinant of the map y i— > X(t,y). 
The problem is that in general there is a (possibly, very short) time T* > 0, at 
which the flow X t ceases to be one-to-one. Points x £ M. d at which this happens 
are caustic points and T* is called the caustic onset time. More precisely, let 

% = {x£ R d : there is y £ R d such that x = X(t, y) and J t (y) = 0}, 

then the caustic set is defined by e ta :— {(x, t) : x £ 'tat} an( i the caustic onset time 
is 

T* := inf{t £ R+ : C t ^ 0}. 
Fort > T* the solution of (1.3), obtained by the method of characteristics, typically 
becomes multi-valued due to the possibility of crossing trajectories, see Fig. [T] On 




Figure 1. Classical trajectories for initial data WSq(x) 
— tanh(5a; — |) 



the other hand, weak solutions to ( 1.3 ), which can be uniquely defined (for example, 
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by invoking the Lax-Olejnik formula) are not smooth in general and thus plugging 
fL2| into ([□} is no longer justified. From the physical point of view T* marks 



the generation of new frequencies within tjj e not captured by the simple one phase 
WKB ansatz ( |1.2| . Indeed, it is well known that for t > T* one generically requires 
a multi-phase WKB ansatz to correctly describe the asymptotic behavior of ip £ , see 
Section [5] for more details. 

1.2. Bohmian trajectories and Bohmian measures. It is well known, that to 
any wave function -0 e £ H 1 (R d ) one can associate two basic observable densities. 
Namely, the position and the current-density, defined by 

(1.9) P e (t,x) = \ip e (t,x)\ 2 , J £ (t,x) = eIm(^(i,a;)V^ e (t,a;)), 
which satisfy the conservation law 

d t p £ + div x J £ = 0. 

These two quantities play an important role in Bohmian mechanics developed in [5J 
|S] (see also [T7] for a general introduction) . In this theory, one defines e-dependent 
trajectories X s : y n- X £ (t,y), via the following differential equation 

X £ (t, y) = u e (i, X £ (t, y)), X £ (0, y) = y £ R d , 

where the initial data y £ R d is assumed to be distributed according to the mea- 
sure Pq = IV'ol 2 •= L+(R d ) and u e denotes the quantum mechanical velocity field, 
(formally) defined by 

^ F , s J e (t,x) fVlJj S (t,x) 

(1.10) u e (t,x):= — y—(=elm' 



Note that the kinetic energy of %p e can be written in terms of p e and u e as 

E k in(t) ■= % [ \V^(t,x)\ 2 dx 



p E (t,x)\u 6 (t,x)\ 2 dx + — / \S/yfff{t,x)\ 2 dx, 



(1-11) , 

_ 2 

which allows to define u e £ L 2 (R d ; p £ dx). Although u e is not necessarily continuous, 
it was rigorously proved in [SJ 135] that, for all t £ K, x 1— > X e (t, ■) is well-defined 
Pq — a.e.. In addition, one finds (under mild assumptions on the potential V) that, 
for all times t £ M + , the position density p E (t,x) is given by the push-forward of 
the initial density Pq(x) under the mapping X t e , i.e., for any non- negative Borel 
function a : R d — >• [0, +00] it holds: 

(1.12) / a(x)p £ (t,x)dx= [ <j{XC{t,y))pl{y)dy. 



In |32j , a phase space description of Bohmian mechanics was rigorously introduced 
through the definition of a class of positive measures, called Bohmian measures, 
(3 s £ M + {R d x R d ), associated to p e and J £ . (Here and in the following, M + 
denotes the set of non-negative Borel measures.) Indeed, for a given scale e > 
and any ip £ £ H 1 (R d ), one defines the associated Bohmian measure f3 £ by 

P e {t,x,p) := p s (t,x)S(p-u e {t,x)), 



where u 6 is defined by (1.10) and S(p — •) denotes the d-dimcnsional delta distri- 
bution with respect to the momentum variable p £ R d . Note that even though 
u £ is not well defined at points where p 6 (t,x) — 0, the Bohmian measure f3 e is. 
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In addition, the zeroth and first moment of f3 £ with respect to p £ M. d yield the 
quantum mechanical particle and current densities, i.e., 



P e (t,x) 



/3 £ {t,x,dp), J £ (t,x) 



p/3 £ (t,x,dp). 



It is shown in [32|, Lemma 2.5] that, for all t £ R+, the measure /3 £ (t,x,p) is given 
by the push-forward of the initial measure 

P £ o(y,p) = P £ o(y)6( P -ul(y)), 

under the following e-dcpcndcnt flow: 

f X*(t,y)=P e (t,y), X £ (0,y)=y, 

{ P%t, y) = -VV(X £ (t, y)) - VVg(t, X s (t, „)), P £ (0, y) = u e (y), 
where Vg(t,x), denotes the Bohm potential: 



(1.13) 



(1.14) 



VS := 



More precisely, for any non-negative Borel function 99 : x — > [0, +00] it holds 



(1.15) 



ip(x,p)/3 £ (t, dx, dp) 



<p(X £ (t,y),P e (t,y))p £ (y) dy. 



Note that (1.131 is the characteristic flow of the following perturbed Burgers' type 
equation 



d t u £ + (u £ • V) u £ + W = W§(t,x), u £ \ t=0 



'o- 



which allows us to identify X £ (t, y) = P £ (t, y) = u £ (t, X £ (t, y)). On the other hand, 
for wave functions ip £ given in WKB form (1.2), we have J £ = \a £ \ 2 WS in which 
case the velocity field is simply given by 

u £ (t,x) = ^\=VS(t,x). 
p £ {t,x) 



One can therefore regard ( 1.13 ) as a nonlinear perturbation of the classical equations 



of motion (1.6). One consequently expects the Bohmian trajectories (X e ,P £ ) to 
converge to the corresponding classical (X,P), in the limit e —> 0+. We shall 
prove that, at least before caustic onset, this convergence indeed holds true (in 
a sense to be made precise, see Theorem |3.1[ ). After caustic onset, however, the 
situation in general is much more complicated in view of Fig. [T] Indeed, we shall 
show that in general one cannot expect the Bohmian trajectories to converge to the 



(multi-valued) classical flow, see Theorem 5.4 



The rest of the paper is organized as follows: In Section [2] we describe some 
general properties of Bohmian dynamics and of the Young measures associated to 
the Bohmian trajectories. These properties will be used in Section [3] to prove that 
the Bohmian trajectories converge to the classical ones before caustic onset. In Sec- 
tion [4] we prove a general result about Bohmian measures associated to multi-phase 
WKB states. This result is then used in Section[5]to show that, even in the free case 
(where V(x) = 0), the Bohmian measure may differ from the Wigner measure, and 
that in general the Bohmian trajectories do not converge to the Hamiltonian ones 
after caustics. Finally, in Section [6] we present a numerical simulations of Bohmian 
trajectories in the regime < e <C 1. 
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2. Mathematical preliminaries 

2.1. Basic a-priori estimates and existence of a limiting measure. From 
now on the potential V will satisfy the following assumptions. 

Assumption 2.1. The potential V € C°°(]R d ;]R) is assumed to be bounded below 
and sub-quadratic, i.e., 

d%V G L°°(R d ) , Vfe e N d such that \k\ > 2. 

Since U is bounded below, without loss of generality we can assume V(x) ^ 0. 
Assumption 2.1 is (by far) sufficient to guarantee the existence of a unique strong 
solution i\) e € C(R d ; L 2 (R d )) to ( |1.1[ ), satisfying two basic conservation laws of 
quantum mechanics. Namely, conservation of the total mass 



(2.1) M £ (t):= I \ip e (t,x)\ z dx = M e (0), 
and the total energy 

(2.2) E £ (t) := %- I \\7^ £ (t,x)\ 2 dx+ I V{x)\ifj £ (t, x)\ 2 dx = E £ (0). 



A direct consequence of these conservation laws is the following result to be used 
later on. 



Lemma 2.1. Let V satisfy Assumption 2.1 and ipQ £ H 1 (M. d ). Then, it holds: 
rT 

P £ (t,y)\ 2 p £ (y)dydt^TE £ (0), VT e E+. 



o 

Proof. Let us recall that p £ (t,x) is the push forward of pg under the mapping X t e , 
i.e., identity (1.12 1 holds true for all t £ M.+ . Using this identity with cr(-) = 
\P £ (t, -)| 2 and recalling that P £ (t,y) = X £ (t,y) = u £ {t, X £ (t,y)), we find 

T I \P £ (t,y)\ 2 p £ (y)dydt = f [ \u £ (t, X £ (t,y))\ 2 p £ Q (y) dy dt 

JR d JO JR d 

T 

2 „e/ 



\u £ {t,y)\ 2 p £ (t,y) dydt. 

JR d 



In view of energy conservation, the last term on the right hand side is bounded by 
\u £ (t,y)\ 2 p £ (t,y)dydt^ [ E £ (t) dt = TE £ (0), 



10 JR d 

as desired. □ 



In addition, to Assumption |2.1| we require the following basic properties for the 
initial datum i/jq. 



Assumption 2.2. The initial data of (1.1 1 satisfy M e (0) = HV'ollia = 1, and there 
exists Co > such that 

sup E £ {Q) < C . 

0<e<l 

Remark 2.2. The normalization HV'ollia = 1 is imposed for the sake of mathe- 
matical convenience. ^From a physical point of view, it is required for the usual 
probabilistic interpretation of quantum mechanics in which p £ — \i> £ \ 2 denotes the 
probability measure of finding the particle within a certain spatial region f2 C M. d . 
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Assumption 2.2 together with conservation of mass and energy and the fact that 
V{x) > 0, implies that for all t e M+: 



(2.3) 



sup (\\*p e (t)\\ L 2 + \\eVi; E (t)\\ L 2) < +oo. 



In other words, ip £ (t) is e -oscillatory and we are in the framework of [32) . Indeed, it 
was shown in j3H Lemma 3.1] that (2.3) implies the existence of a limiting measure 
/3(t) £ A4 + (R d x Mz) such that, up to extraction of a subsequence, it holds: 



(2.4) p 
and we also have 



e-s-Oa 



(3, in L°°(IR t ;X + (IR^ x IE^))weak- 



(2.5) 



p e (t,x] 



£->0 4 



P(t,x,dp), J E (t,x) 



pfi(t,x,dp), 



where the limits have to be understood in L co (M. t ; A4 + (M. d )) weak— *. 

2.2. Young measures of Bohmian trajectories. The limiting Bohmian mea- 
sure /3 is intrinsically connected to the Young measure (or parametrized measure) of 
the Bohmian dynamics. To this end, we first note that <& e (t, y) = (X s (t, y),P e (t, y)) 
is measurable in t, y and thus, there exists an associated Young measure 



T 



t,y 



M+(R d x Rp) : (t,y) ^ T t , y (dx,dp), 



which is defined through the following limit (see [4, 25, 35 ): for any test function 



lim 

£->0 



a(t,y,<S> e (t,y))dydt 



a(t,y,x,p)T t , y {dx,dp) dydt. 



Having in mind (1.15), if we assume in addition that 

Po ^ Po , strongly in h\ (R d ) 
we easily get the following identity: 



(2.6) 



P{t,x,p) 



^t,y{x,p)p (y)dy. 



Here, /3 is the limiting Bohmian measure obtained in ( 2.4 ) for a specific subsequence. 
The relation (2.6) has already been observed in |32j and can be used to infer the 
following a-priori estimate on T t>y . 

Lemma 2.3. Let Assumptions \2.1\ and \2.2\ hold, and assume in addition that 

Pq £ ^^~ po strongly in L^_(M. d ). Then, for any T G R + . there exists a C — C(T) > 
such that 

\p\ 2 Po(y)T t , v ( d ^^P)dydt^C(T). 



Proof. Using ( |2.6| we see that 

\p\ 2 Pa{y)^t, y {dx,dp)dy 



\p\ 2 (3(t,dx,dp). 



Now we recall that, by definition, 

p e (t,x,p)=p £ (t,x)5(p-u E (t,x)) 

and hence 



\p\'P s (t,dx,dp) 



p e (t,x)\u'(t,x)\ 2 dx < 2E% a (t) < C(T), 
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in view of ( |1.11[ ) and energy conservation. This uniform (in e) bound together with 
Fatou's lemma implies 



\pff3(t,dx,dp)^C(T), 
and the assertion is proved. □ 



Lemma 2.3 together with Lemma 2.1 will be used to prove the following impor- 
tant property for the zeroth moment of T t)y . 

Proposition 2.4. Let Assumptions\2J\ and\2J^ hold, and assume in addition that 



Pq — > po strongly in L^_(R ). Denote 

v t , y {x) := I Y t , y (x,dp). 



Then Vt, y G A4 + (BL X ) solves, a.e. with respect to the measure po(y), the following 



transport equation 



d t v t ,y + div x (^J pT t , y (x, dp)^j = 0, v t =a, y {x) = S(x - y), 

in the sense of distributions on D'(M. t x M. x ). 

This transport equation will play a crucial role in the convergence proof of 
Bohmian trajectories before caustic onset. 

Proof. As a first, preparatory step we shall prove that, for all test functions C G 
C (R t x Ry), a G C {R d x ): 



6-S-O 



+ Jo 



lim / / P £ (t 7 y)at,y)a(X s (t,y))f (y)dydt = 



C(t,y) // po-(x)T t , y (dx,p) p (y)dydt, 



(2.7) 



To this end, let K > and X i< G C c °°(IR d ) be such that and xk(p) = 1 for \p\ < K, 
and xk(p) = for \p\ > K + 1. Then, by writing P e = Xk{P £ ) + (1 - Xk(P 6 )) we 
can decompose 

C(t, y)o-(X%t, y))P s (t, y)pl(y) dydt = I[' K + % K . 

Because of the strong convergence of pg, the first term on the right hand side has 
the following limit: 

Ii KS —f C(t,y) a(x)xK(p)Tt, y (dx,dp) p (y)dydt, 



On the other hand, by having in mind the result of Lemma |2.1| the second term 
on the right hand side can be estimated by 



I^Kcf/ \P^(t,yM(y)dy, dt 

JO J\P'\^K 

/ / \P e (t,y)\ 2 pl(y)dy,dt<^-E s (0). 

& JO JR* K 



In view of Lemma 
validity of p7|. 



2.3 



we can let K — > +oo, which yields A ^>°° Q and the 
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With (2.7 1 in hand, we shall now show that vt, y indeed obeys the transport 
equation given above. Let £, <p G C£°(R d ), a G C^°[0, oo), be smooth compactly 
supported test functions. Then by (2.7) we get 

(d t a(t)<p{x) + a(t)p ■ V x ip(x)C(y))T yit (x,dp)p (y) dydt 

2d V / 

(dMt)<p(X e (t,y)) + a(t)P e (t,y) ■VMX e (t,y))C(y))p (dy)dt. 



= lim 



+ Jo 

Recalling that P £ (t,y) = X e (t,y), which implies that 

P £ (t, y) ■ WMX e (t, y)) = j t v{X s {t, y)), 

we obtain 

\ihoM) s> X (7 . // ) ! + a ( / i Pt(t, y) • V , vl A" : (/.,/) K ( i ) ,, u (,/.,/ ) ,// 
i( V,.rT(/) y - ( A-(/.//ii- 'Ti'i'' 

/O ^R d 

a(0)^(X e (0,2/))C(y)po(?/)d2/ 



JR d 



(d t a(tMX £ (t, y)) + a(t)j t p(X s (t, y))C(y)) p (y) dy 



dt 



•^(o) l p(y)C(y)po{y)dy. 

where in going from the second to the third we have itegrated by parts with respect 
to time, and from the third to the forth line we have used that X £ (Q,y) = y 
by definition. The obtained expression in the last line is nothing but the initial 
condition, since 



// 



( P(y)C(y)po(y)dy = // ip{x)T , y (x,dp)C{y)p {y)dy, 

J jR2d 

is equivalent to saying that 

v , y (x)= / T 0i y(x,dp) = 5(x - y), p (dy) - a.e. 

JR d 

□ 

Having collected all necessary properties of Y t , y we shall prove the convergence 
of Bohmian trajectories (before caustic onset) in the next section. 

Remark 2.5. For completeness, we want to mention that T t . y is indeed a proba- 
bility measure on x for a.e. y, t, provided the sequence {V' e }o<e^i is compact 
at infinity (tight), i.e., 



lim limsup / \tfj £ (t,x)\ 2 dx = 0. 



Indeed if the latter holds true, it was shown in [351 Lemma 3.2] that 



limM e (t)= lim // (3 £ (t,dx,dp) = // P(t,dx,dp), 

£ ^°+ £->0+ J J M 2d J J Kid 

and having in mind our normalization M £ {t) = 1, we conclude 

1= // P&dx,dp)= p (y)T ttV (dx,dp)dy, 

J JR 2d JJjR 2d xR d 



in view of (2.6). Define 

a t,y : = / / Tt,y(dx,dp) < 1. 
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Then, since J Rd po{dy) = 1, we conclude a Vt t — 1 a.e.. However, we shall not use 
this property in the following. 

3. Convergence of Bohmian trajectories before caustic onset 

So far we have not specified the initial data ipo to be of WKB form. By doing so, 
we can state the first main result of our work (recall the definition of sub-quadratic, 



given in Assumption 2.1 1. 

Theorem 3.1. Let Assumptions \2A\ hold, and let i/jq be given in WKB form 
(3.1) r (x) = a (xy s °^/ £ , 

with amplitude oq G 5(IR <i ;C) and sub-quadratic phase So G C°°(]R d ;IR). Then, 
there exists a caustic onset time T* > such that: 

(i) For all compact time-intervals It C [0, T*), the Bohmian measure (3 e associ- 
ated to p e , J £ satisfies 

(3 £ p(t,x)S(p-WS(t,x)), in L°°(I t ;M + (Ri x R d ))weak-*, 



where p G C°°{I t ;S{R d )) and S G C°°(I t X R d ) solve the WKB system |L5| , ( [L31 ). 
(ii) The corresponding Bohmian trajectories satisfy 

X* e ^?X, P^^ + P 

locally in measure on {It x supppo} C K t x R d , where po = \ao\ 2 . More precisely, 
for every S > and every Borel set Q C {/ ( x supp po} with finite Lebesgue measure 
5? d+1 , it holds 

lim J? d+1 ({(t,y) G Q : |(X e (t, y), P%t, y)) - (X(t, y), P{t, y)) \ > 8}) = 0. 

Assertion (i) above was already proved in [32], but since the obtained form 
of the limiting measure will be used to show (ii), we shall recall the proof of (i) 
for the sake of completeness. Assertion (ii) shows, that before caustic onset, the 
Bohmian trajectories converge locally in measure to the corresponding classical 
Sow. Clearly, if ao(x) > for all x £ K d , and thus supppo = we obtain local in 
measure convergence of the Bohmian trajectories on all of J t xKf. After selecting an 
appropriate sub-sequence {s n }nefi this also implies (see, e.g., [7]) almost everywhere 
convergence on any finite subset of I t xM. d . Moreover, since, by definition, X s = P e , 
the convergence in measure of P e to P combined with the L 2 bound from Lemma 
2.1 implies that, for .5f d -a.e. y, the curves X £ (-,y) converge uniformly to X(-,y) 



on the time interval It- 



Proof of Theorem \3.1\ We first note that (3.1 ) implies 

_2 



E £ (0) = l f \a \ 2 \VS\ 2 dx+ S - ( \Va \ 2 dx+ f V(x)\a \ 2 dx. 

* JR d 1 JR d JR d 



Since a G S(R d ), we see that Assumption 2.2 is satisfied and thus all the results 
established in Section [2] apply. In particular, we have the existence of a limiting 
Bohmian measure /3 G L°°(Rt;M + {Ri x R d )) weak—*. In order to prove Assertion 
(i) we need to show that before caustic onset, this limiting measure is given by a 
mono-kinetic phase space distribution, i.e., 

(3.2) (3(t,x,p) = p(t,x)6(p-VS(t,x)). 

In [3 2) sufficient conditions for (3 being mono-kinetic have been derived. In partic- 



ular, it is proved in there that (3.2 1 holds as soon as one has strong L 1 convergence 



of p E and J e in the limit e — > 0+. To show that this is indeed the case, we shall 
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rely on the so-called modified WKB approximation introduced in [23] and further 
developed in [10]: Define a complex- valued amplitude a £ by setting 

(3.3) a £ (t, x) = i> £ (t, x)e~ lS{t ' x)/e , 

where ip £ solves (1.1) and S is a smooth solution of the Hamilton- Jacobi equation 

(1.3) . Next, we recall that the results of [TU] (see also [TT]) ensure that under 
our assumptions there is a time T* > 0, independent of x <G R d , such that, for 
all compact subsets It C [0,T*), the Hamiltonian flow (1.6) is well-defined, and 
there exists a unique (sub-quadratic) phase function S £ C°°(I t x M d ), given by 
(1.7). Consequently, this also ensures the existence of a smooth amplitude a G 
C^iltlSiW 1 )) given by ([□§. 

With this result in hand, a straightforward computation shows that a £ , defined 
in (3.3), solves 

(3.4) d t a £ + Va e • VS + °^AS = «|Aa e , a £ (0, x) = a Q (x). 



This equation can be considered as a perturbation of (1.4). Indeed, if we denote 
the difference by w £ := a £ — a, then w £ satisfies 

d t w £ + \7w £ • VS + ^- AS = 4 Aa '' U ' e ( ' ^ = °' 

where the source term on the right hand side is formally of order 0(e). Invoking 
energy estimates, one can prove (see [101 Proposition 3.1]) that for any time-interval 
It C [0,T*), there exists a unique solution a £ G C(I t ; H s (R d )) of fl3~I| ), and that 

\\w e \\L°°(I t ;H*(BL-')) = \\ a£ - a\\L°°(I t ;H*(R*)) = 0(s), Vs ^ 0. 

Writing the mass and current densities as 



p £ = \,jj £ \ 2 = \ a £ \ 2 , J £ = elm^V^) = |a e | 2 V5 + eIm(a E V 



ar 



and using the fact that H s {W i ) L°°(M. d ) for s > d/2, this consequently implies 
p £ p, in L°° (I t ; L 1 (R d )) strongly, 

and 

J £ ^ pu, in L°° (It ; L\ oc (R d ) d ) strongly, 

where p = \a\ 2 and u = V5 are smooth solutions of the WKB system: 

d t p + div x (pu) = 0, p(0,x) = \a (x)\ 2 , 

d t u + u-Vu + VV(x) = 0, u(0,x) =VS (x). 

In particular, we infer that P(t,y) = V S(t, X(t,y)) — u(t, X(t,y)) and, in view of 
(1.8), we also have that the density p = \a\ 2 is given by 

3-5 p(t,x) = ———— te[0,T). 

J t (Y(t,x)) 

The strong convergence of p £ , J £ together with [321 Theorem 3.6] then directly imply 
that the limiting measure f3 is given by (3.2) and thus Assertion (i) is proved. 

In order to prove (ii) we first note that for every fixed t £ [0,T*), the limiting 
measure j3(t) is carried by the set 

<S t = {(x,p) G R 2d :p = u(t,x)}. 



The identity (2.6) then implies that a.e. in y the measure T t . y is also carried by 
the same set and we consequently infer 

Y t ,y(x,P) = Ht,y(x)6(j3 ~ U(t, X)), 

where p t , y is the Young measure associated to X £ (t,y). 
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By taking the zeroth moment of Tt, y with respect to p £ M. d we realize that 
indeed pt, y — Vt,y, with v t , v defined in Proposition 2.4 We thus find that p t ,y 
solves, in the sense of distributions: 

dtVt, y + divz (ufi t , v ) = 0, fH=o,y(x) =6(x-y), 

a.e. with respect to the measure po{y). In other words, fj,t, y {x) solves the same 
transport equation as the limiting density p(t, x) does. In view of ( |3.5[ ), we therefore 
conclude that, before caustic onset, /x tjl/ is given by 

Ht, y (x) = jtyn X )) 5 ( Y ^ x ) ~ y)' P° ~ a - e - 

Multiplying by a test function ip £ Co(K-f x and performing the change of 
variable x = Y(t,x), we consequently find 

(Vt,y, ( P)=[[ T rvT. — ^§(Y(t,x) -y)<p(x,y)dxdy= [ <p(X(t,y),y)dy, 
JJ R2d J t {Y{t,x)) J Rd 

and thus we can also express fit.y — S(x — X(t, y)). In summary we obtain that 

T t Jx, P ) = S(x - X(t,y))S(p - u(t,X(t,y))). 

a.e. on supppo Q R • In other words, the Young measure T 4)?/ is supported in 
a single point (on phase space). By a well known result in measure theory, cf. 
[25| Proposition 1], this is equivalent to the local in-measure convergence of the 
associated family of trajectories X £ , P £ and we are done. □ 

The proof in particular shows, that, at least before caustic onset, the Young 
measure T t . y is independent of the choice of pg, even though the Bohmian flow X s 
is not. 



Remark 3.2. It is certainly possible to obtain Theorem |3.1| under weaker regularity 
assumption on V, Oq, and So, which are imposed here only for the sake of simplicity. 
The assumption of V and So being sub-quadratic, however, can not be relaxed, if 
one wants to guarantee the existence of a non-zero caustic onset time T* > 
uniformly in x £ M d , see, e.g., [10] for a counter-example. Explicit examples of 
initial phases So, for which T* = +oo (i.e., no caustic) are easily found in the case 
V(x) = 0. Namely, either plane waves: Sq(x) = k ■ x, where k £ M. d is a given 
wave vector, or Sq(x) — — \x\ 2 , yielding a rarefaction wave for t £ K + , see [T^] . 
In these situations, we obtain in-measure convergence of the Bohmian trajectories 
(X s , P £ ), and consequently also uniform convergence of X £ , locally on every Borcl 
set C {K t x supppo} with finite Lebesgue measure. 



4. Superposition of WKB states and Bohmian measures 

4.1. Bohmian measure for multi-phase WKB states. In view of Fig. [T] we 

expect that for \t\ > T* , i.e., after caustic onset, the correct asymptotic description 
of ip £ is given by a superposition of WKB states, also known as multi-phase ansatz. 
In order to gain more insight in situations where this is indeed the case we shall, 
as a first step, study the classical limit of the corresponding Bohmian measure. To 
this end, let £1 C R t x be some open set and consider ip £ to be given in the 
following form: 

N 

(4.1) r(t,x) = Y,bj(t,xy Si(t ' x)/£ +r e (t,x), 

3=1 
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where bj <E C°°(f2;C) are some smooth amplitudes and the real- valued phases 
Sj € C°°(n;R) locally solve 

1, 



(4.2) 



d tSj + - \VSj\ J + V(x) = for all j = 1, 



,N, 



In addition, r £ denotes a possible remainder term (the assumptions on which will 
be made precise in the theorem below). 

Remark 4.1. As we shall see Section|5] the multi-phase WKB form (4.1 1 can be 
rigorously established, locally on every connected component of (R< x R*) \ ^, i.e., 
locally away from caustics. 

The second main results of this work establishes an explicit formula for the 
limiting Bohmian measure /3 associated to a wave function of the form (4.1 ). More 
precisely we prove the following: 



Theorem 4.2. Let ip e be as in pi) , with bj € C°°(ft;C), Sj £ C°°(ft;R), for all 
j = 1, . . . , N, where ft C [0, T\ x R d denotes some open set. Assume, in addition, 

(4.3) VS^VSk for all j^k £{!,..., N}, 

and that the remainder r £ (t,x) satisfies 

(4-4) \\r s \\ L?oc{n) =o(l), \\eVr e \\ L?oc{n) = o(l) ase^0+. 

Then 

(3 E f3(t,x,p), in L°°([0,T};M + (M. d x x R^))weak-*, 
where j3 is given by 



P(t,x,p) 

with 6 = 

(4.5) T(t,x, 



6 N ) E T N , ano 

N 



r , , ,, Ef, k ^Sj(t,x)T jik (t,x,e) 
L(t,x,0)o[p — 



T(t,x,8) 



d6. 



E 6 ^ 

3=1 



x e 



it). 



■3,k 



(t,x,9) :=Re {b 3 b k e i{e ^ ek) ^ 



The above formula for (3 generalizes equation (6.6) given in [35] and states that 
j3 in general is a diffuse measure in the momentum variable p € R d , unless all but 
one bj — 0. Note that, in the case where N — 1, f3 simplifies to a mono-kinetic 
phase space measure, i.e., 

f3(t,x,p) = \b(t,x)\ 2 6(p-VS(t,x)). 



We already know from Assertion (i) of Theorem 3.1 that this holds for \t\ < T* 
i.e., before caustic onset. 



Proof. By our assumptions, it is easy to check that p 6 — \ip 
J 6 = eIm(V^(t,x)VV' e (i,a;)) = J 6 + r 2 , e , where 



,£ — U/,e|2 _ X£ 



r\ >e and 



N 



N 



p e := b~b k e^ s '- Sk),e , J e := ^ VSj Re (b~b k e l(s ^ Sk)/e 



and 



lki, s |Ui oo (n) = o(l), lk 2 ,e||i,i oc (n) = o(l). 
In order to derive the classical limit as e — > 0+ of the Bohmian measure /3 s , we 
need to compute the limit of expressions of the following form 

J e {t,x)- 



(4.6) 



a(t,x)p (t,x) V {t,-^—j 



dx dt, 
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where ip, a <E C^°([0, T] x M. d ; K) are smooth and compactly supported. To this end, 
we first note that, because tp is smooth and compactly supported, the map 

R+ x R d 3 (s,v) i ^ stp(t, 

is Lipschitz (uniformly with respect to t), which implies 

^c(\\r he \\ Ll (n) + ||r 2 , e |Ui (o) 



) =0(1). 



In particular, to compute the limit as e — > 0+ of the expression in (4.6) it suffices 
to consider 



(4.7) 



/ J s (t 

y(t,x)p e (t,x)ip[t, — - 1 



dx dt. 



J s {t,xy 

,x). 

We now use the following result, whose proof is postponed to the end. 

Lemma 4.3. There exists a set S C f2 of J2f d+1 -measure zero such that, for all 
j,k,£€ {1,...,N} with k^i, 

Sj(t,x) - S k {t,x) 



Sj(t,x) - S t (t,x) 



for all {t, i)6fi\E. 



Using this lemma, we deduce that for jf d+1 — a.e. (t,x), the frequencies 
Si(t,x) - S k (t,x) 



k = 2, 



,N, 



are all rationally independent, which implies that the "trajectories" 

' S2 — >Si * 



e 1 ^ 



( ( l 
I cos I - 



( Sn — Si 
, cos^ 



and 



£ 1 — ^ sm 



fS 2 -S 1 



, sm 



Sn — Si 



are both dense on the (N — l)-dimensional torus T^ -1 . By standard results on 
two-scale convergence (see for instance [I]), we consequently obtain that for any 
continuous and compactly supported test function $ : 51 x C^ -1 — > K, 



J $(t,x,e^- s ^,...,e^- s ^/ £ ) dxdt 



Moreover, we observe that for any j, k we can write 

Sj — Sk Sj — Si Si — Sk 



(t, x, e ldl , . . . , e^"- 1 ) de x ... d6 N -i dx dt. 



Hence the expression in (4.7) converges to 



N 

r(t,x) / V b b k e l 



if\ t, 



Ef*=iV5 i Ite(6 i 6 k e < 



E^,fe=i bjb k eWi-*- $ >-0 



d9i . . . d0N-i dx dt, 



where by convention 9q = 0. Finally, let us observe that one can also rewrite the 
obtained expression in a more symmetric form by performing the change of variables 
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6j-i -H- Oj — 6±, and it is immediate to check that under this transformation the 
above expression is equal to 

cril.r) I V It, X, 9) <p [t, °' T u X Z dx dt, 



where 6 = (6*i, . . . , 9 if), and T and Tj tk are defined in (4.5). By the arbitrariness of 
ip and a, this proves the desired result. □ 

We are now left with the proof of Lemma [473] 



Proof of Lemma 4-3 The set £ can be described as 

U U m% 

j.k.l, k^l m^neZ 

where 

Sj!Ci-={(t,x)efl : m[S j (t,x)-S k (t,x)} + n[S j (t,x)-S i (t,x)]=0}. 

We now claim that each Sj 1 ^ is a smooth hypersurface in J7, which implies in 
particular that SJ^" (and so also S) has measure zero. To prove that this is indeed 
the case, it suffices to check, in view of the implicit function theorem, that the 
gradient of the function 

(t, x) i-> m[Sj(t, x) — Sk(t, x)] + n[Sj(t, x) — Si(t, x)] 

is nowhere zero. Assume by contradiction that this is not the case, i.e., there exists 
a point (t, x) € Q where 

(m + n)d t Sj(t,x) = md t Sk(t,x) + nd t Si(t,x), 

(m + n)VSj(t, x) = mVSk(t, x) + nVS((t, x). 
By (4.2), the first equation above becomes 

(m + n)\VSj(t,x)\ 2 = m\VS k {t,x)\ 2 +n\VS l {t,x)\ 2 , 
which combined with the second equation gives 



[m + n) 



VS k (t,x) + — —VS t (t,x) 



m\VS k (t,x)\ 2 +n\VS e (t,x)f 



m + n m + n 

By strict convexity of | • | 2 , the above relation is possible if and only if WS k (t, x) = 
VS?(t,x), which contradicts (4.3) and concludes the proof. □ 

4.2. Comparison to Wigner measures. An important consequence of Theo- 
rem [472] concerns the connection between the limiting Bohmian measure j3 and the 
Wigner measure w € M + (M.^, x K~) associated to ip £ . To this end, let us first recall 
the definition of the e-scaled Wigner transform w £ given in [21 [201 (SHI : 

W *(t, x, p) := JL_ £ r (t, (t, x + e - v ) d V . 

Provided il> £ (t) is uniformly bounded in L 2 with respect to e, it is well known that, 
cf. [2ni[3nj there exists a limit w(t,x,p) such that 

w £ w, in L°°(R t ;M + {Ri x R^))weak-*. 

In addition, one finds w(t) 6 M + (M.^xR^), usually called Wigner measure (or semi- 
classical defect measure). The latter is known to give the possibility to compute 
the classical limit of all expectation values of physical observables via 

lim {if) e (t), Op s (a)ip e (t)) L 2 ^ = a(x,p)w(t,x,p)dxdp, 
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where the Op e (a) is a self-adjoint operator obtained by Weyl-quantization of the 
corresponding classical symbol a e S(R d x R^), see [20l[36] for a precise definition. 
In addition, if tp s (i) is e-oscillatory, i.e., satisfies (2.3), we also have that the zeroth 
and first p-moment of w yield the classical limit of p s and J e , i.e., 

p £ (t,x) — ? / w(t,x,dp), J £ (t,x) — ? / pw(t,x,dp), 

JR d JR d 

where the limits have to be understood in L°°(R t ; A4 + (M. d )) weak— *. Note that 
this is indeed analogous to (2.5). 



For a given superposition of WKB states such as (4.1), the associated Wigner 
measure has been computed in [30] (see also [36 ): under the same assumption on 
the phases, i.e., VSj ^ VS^ for all j ^= k, one explicitly finds 

N 

(4.8) w(t,x,p) = \b 3 {t^)\ 2 5( P -VS 3 {t,x)). 

3 = 1 

^From this explicit formula we immediately conclude the following important 
corollary. 



Corollary 4.4. Let bj ^ 0. Then, under the same assumptions as in Theorem 
we have that, in the sense of measures, j3 = w if and only if N = I . 



Proof. For bj ^ and iV > 1 we see from Theorem 4.2 that (3 is a diffuse measure 
in the momentum variable p £ M. d , and thus (3 ^ w in view of (4.8). On the other 



hand, if N = I then, both w and f5 simplify to the same mono-kinetic phase space 
distribution. □ 



In view of Assertion (i) of Theorem 3.1 we conclude that before caustic onset, 
the classical limit of all physical observables can be computed by taking moments 
of the limiting Bohmian measure, since in fact f3 = w for \t\ < T* . After caustic 
onset, however, this is in general no longer the case (see Section j5J). 

Still, we do know (by weak compactness arguments) that the zeroth and first 
moments w.r.t. p £ M. d of /3 and w are the same for all times i € R. For complete- 
ness, we check this explicitly in the case of multi-phase WKB states: using the fact 
that 



we compute 



/ cos (6>) d9 = [ sin (0) d9 = 0, 
Jt Jt 

N 

P{t,x,dp)= V bj(t,x)b k (t,x) / e 1 ^- ^ d6L-.de 

N 

= J2\b J (t,x)\ 2 = / w(t,x,dp). 



Moreover 



/ P (3{t,x,dp) = V VSj(t,x) [ Re(b j (t,x)b k (t,x)e l{ ^- 8k) ) dQ 1 ...dQ 

N 

= J2VS j (t,x)\b j (t,x)\ 2 = / pw(t,x,dp). 

~1 JR d 



In other words, in the case of multi-phase WKB states, the difference between w 
and (3 can only manifest itself in p-moments of order two or higher. 
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5. A COMPLETE DESCRIPTION IN THE FREE CASE AND POSSIBLE EXTENSIONS 

In this section we shall give a (fairly) complete description of the classical limit 
of Bohmian dynamics in the case of the free Schrddinger equation corresponding to 

V = 0. The proof will rely on classical stationary phase techniques. For the case 

V =/= decisively more complicated methods based on Fourier integral operators 



have to be employed, as will be discussed in Section 5.3 



5.1. Multi-phase WKB for vanishing potential. Consider the free Schrodinger 
equation with WKB initial data: 

(5.1) ied t ^ + yAV e = 0, V e |t=o = ao(x)e lS °^/ e , 

In this case, we find the free Hamilton-Jacobi equation, which is obviously given by 

(5.2) dS+ ±\VS\ 2 =0, S\ t=Q = So, 



and the corresponding classical Hamiltonian equations (1.6 1 simplify to 

f x(t,y) = P(t,y), x(P, v ) = y, 

< "" >1 \P(t,y)=0, P(0,y) = VS o (y). 

This implies that, for all t G R+, P(t,y) = V5 (y) and 
(5.4) X(t,y)=y + tVS (y). 

Consequently, the caustic set is given by ^f rcc = {(x,t) : x € ^ t free } where for 
x = X(t, y) we set: 



{x G R d : 3 y e R d satisfying and det(Id + tV 2 S (y)) = 0}. 



In particular, we see that in the free case, the caustic onset time T* > is solely 
determined by the (sub-quadratic) initial phase Sg(y). In order to proceed we need 
to slightly strengthen our assumption on the initial phase So- 

Assumption 5.1. The initial phase So G C°°(M. d ; K.) is assumed to be sub-quadratic 
and 

lim ^#1 = 0. 

lv|->oo \y\ 

In other words we need that Sq grows strictly less than quadratically at infinity. 
This is the same assumption as in [5], guaranteeing that the map y i— > X(t, y) is 
proper and onto. 

In the following we shall denote by x t— > y = Y(t, x) the inverse mapping of 



(5.4) . Clearly, for \t\ > T* this inverse will not be unique in general, i.e., for 
each fixed (t,x) € K f x Kf there is N(t,x) G N and corresponding Yj(t,x), with 
j = 1, . . . , N(t, x), satisfying the implicit relation 

(5.5) Y j (t,x)+tVS Q (Y j (t,x))=x. 



Assumption 5.1 guarantees that in each connected component of (M t x Kr) \ % Ge 
there are only finitely many {Yj(t, x)}. (This follows by properness of the charac- 
teristic map and the implicit function theorem, see [5j Lemma 1.1].) In addition, 
in each such connected component N(t, x) — const. Moreover, under the same 
assumptions on So, we already know that the caustic onset time T* is positive, and 
thus there is exactly one connected component Slo of (Rt x M. d ) \ ^ft-ee containing 
{t = 0}. 
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In order to proceed further, we also recall that the solution of (5.1) admits an 
explicit representation in the form of an e-oscillatory integral 

d 

a (y) e ^(*^)/ e dy, 



(5.6) r(t,x) = 

where the phase is given by 



1 



\j2ixiet 



(5.7) 



$(t,x,y) := S (y) 



\x - yf 
2t 



5.1 



It is well known that, for e — > 0+, the representation formula (5.6) can be treated 
by the stationary phase techniques (see, e.g., Theorem 7.7.6. of |24|) and we con- 
sequently obtain the following lemma. 

Lemma 5.1. Let a € <S(R d ;C) and Sq satisfy Assumption 
(t,x) G (M. t x R d ) \ < ^f ree the solution of (5.1) satisfies 

(5.8) ^(t,x) E ^= + a j {t,x)e iwK ^/ A e i ^ t ' x ' Y ^ t ' x ^/ e + r e (t,x), 

3 = 1 



Then, for all 



where $(t,x,y) is given by (5.7), Kj(t,x) € N denotes the Maslov factor, and 



(5.9) 



CLj(t, X) 



a {Y 3 (t,x)) 



|det(Id + tV 2 5 (r,(t,x)))| 1 /2- 
In addition, the remainder r e satisfies 

(5.10) ||r e ||co ( n) = 0(e), ||r e || c i (n) = O(l) as e -> 0+, 

uniformly on compact subsets f2 C (Mt x K d ) \ "^free- 

Remark 5.2. The first remainder estimate ||T e ||c°(f2) = 0( e ) is classical, whereas 
the second one can be obtained by noticing that the operator V commutes with 
the free Schrodinger equation (5.1). Thus, we find that VV> £ satisfies an integral 
representation analogous to (5.6), i.e., 

V^ e {t,x 



/ 1 



J\x-y\ 2 /(2te) V ^(y)cfy. 



V \Jlmet, 

By applying the stationary phase lemma to this oscillatory integral one readily 
infers the estimate HreHc^Mi) = 0(1)- 



Next, we note that, in view of (|5.7| and (|5.5|), we explicitly have 
(5.11) 



$(t,x,Yj(t,x)) = S (Y {t,x)) + j t \x -Y 3 (t,x)f 



S (Y j (t,x)) + ^WSoiYj^x)]' 



On the other hand, since for V(x) — it holds that P(t,y) = VSo(y) (i.e., P is 
cc 
1, 



constant along the characteristics), the solution formula (1.7) yields, for all j 



S 



(5.12) 



1 



•jfax) =S (Yj(t,x))+ -\P(r,y)\ 2 dr\ (t<x) 



= S (Y j (t,x)) + -\VS (Y j (t,x)\ 2 . 



We consequently infer that $(t, x, Yj(t, x)) = Sj(t,x) is a smooth solution of the 
free Hamilton- Jacobi equation (5.2) for all j = 1, . . . ,N(t,x). Obviously, we also 
have that a 3 given by (5.9) solves the corresponding transport equation (1.8) with 
S = Si. 



18 



A. FIGALLI, C. KLEIN, P. MARKOWICH, AND C. SPARBER 



Remark 5.3. An alternative way of showing that $>(t,x,Yj(t,x)) solves the free 
Hamilton- Jacobi equation is to plug (5.11) into (5.2) and use (|5.5|) to implicitly 



differentiate with respect to t and x. A lengthy but straightforward computation 
then yields the desired result. 

For completeness we also recall that the Maslov factor is explicitly given by [21] 

N 3 Kj(t, x) = m~j(t, x) — mj(t, x), 

where m (t, x) € N denotes, respectively, the number of positive or negative eigen- 
values of the matrix Id + tV 2 So(Y 3 (t, x)). Note that Kj can also be written in the 
form 

Kj (t,x) = d — ZmJ (t,x). 
By the implicit function theorem, n(t, x) — const in every connected component of 
(Kt x R^)\^ frGO , see, e.g., 0. 

5.2. WKB analysis of Bohmian dynamics in the free case. From what is 
said above, we infer that in each connected component SI of (Kt xEj)\ ^frcc, the 
solution tp £ admits the approximation (5.8), so Theorem 4.2 can be applied after 
identifying 

bj(t,x) = a^xy™^'* = a j (t,x)e i " ta ^, j = l,.. -,N(t,x) = N a , 

where kq € K and Nq G N are constants depending only on S7. Consequently, we 
obtain the following result. 



Theorem 5.4. Let a Q G <S(R d ;C) and Sq satisfy Assumption \5. l\ Denote by Qq 
the connected component of (Rf x K^) \ ^f ree containing {t = 0}. Then it holds: 

(i) The limiting Bohmian measure satisfies 

/3(t,x,p) = w(t,x,p) = p(t,x)5(jp-u(t,x)), V(t,x) € Sl n , 
and the Bohmian trajectories converge 

X°(t,y) ^°>+ y + tVS Q (y), P%t,y) ^°> + VS (y), 

locally in measure on f2o H {Kt x supppo}- 

(ii) Outside ofQo there are regions fl C (K t X m*) \ ^ frce where (5 ^ w and where 
the Bohmian momentum P E does not converge locally in-measure to the classical 
momentum P. 

(iii) There exist initial data ao(y) and So(y) such that, outside of £1$, there are 
regions Cl C (Kt x K^) \ ^froc in which both X s and P £ = X 6 do not converge to 
the classical flow. 



Note that Assertion (i) is slightly stronger than Theorem 3.1 (i) in the sense 
that Slo is strictly larger than [0, T*) x K^. The proof shows that if |ao| > on all 
of R d , Assertion (ii) holds for any connected component £1 7^ Slo, whose boundary 
intersects the boundary of Sl . 

Proof. We first note that for all (t, x) 6 Qq it holds N(t, x) — 1 and Kj(t, x) = 0. In 
view of the remainder estimates stated in Lemma |5.1| we thus can apply Theorem 
lO with N = 1 to obtain 



P(t,x,p) =p(t,x)S(p- VS(t,x)), 

where p = \a\ 2 . With this in mind, the result on the convergence of the Bohmian 
trajectoriess follows verbatim from the proof of Theorem 3.1 (ii). This proves the 
first assertion. 

In order to prove Assertion (ii) , we first note that that outside of S7o we have (in 
general) more than one branch, i.e., N(t, x) > 1. For instance, assume that |ao| > 
on K d , and let O 7^ S7q be a connected component whose boundary intersects the 
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boundary of Oq. Then it is not difficult to see that Nq ^ 1, as otherwise one could 
show that no caustics can occur on d^loDdU. Next, we recall that in each connected 
component of (R t x M. d ) \ % Tee the phase x, Yj(t, x)) = Sj(t,x) is a smooth 
solution of the Hamilton- Jacobi equation (1.3). By the method of characteristics 
we have that 

V^(t,x,Yj(t,x)) = VSj(t,x) = P{t,Yj(t,x)) = VSoiYj^x)), 

since P(t, y) is constant along characteristics (recall that V(x) — 0). Hence, assum- 
ing by contradiction that VS^ = VS^ for some j ^ k , the above identity together 



with ( |5.4[ ) yields Yj(t,x) = Yfc(i, x), which is impossible by construction. This im- 
plies that in each connected component f2 we can apply Theorem |4.2| to conclude 
that /3 in general is a diffuse measure in p 6 E d , unless all but one of the ctj = 
in fi. In view of (5.9), the latter cannot be the case if |ao| > on 



Corollary 

|4.4| then immediately implies (5 ^ w. On the other hand, since for WKB initial 



data we have that p £ is indeed e- independent, we can apply (2.6) in Vl to infer that 
the Young measure T t . y is diffusive in p (since /3 is). This, however, prohibits the 
convergence of P £ locally in measure, since the latter is equivalent to T t)J/ being 
concentrated in a single point. 

The result in (ii) may still give some hopes for the convergence of X £ to X, since 
the fact that X £ — P e gives more compactness for the curves in the ^-variables. 
However, we shall see that this is not the case. 

Consider indeed the example described in Fig. [ljand Fig. [7] (so d — 1). These 
figures suggest that for i(jq as in (6.5) convergence should not hold. To show this 
rigorously, we begin by observing that p(t, x) > on l f x M. x (this follows from 
the explicit formula for p = \a\ 2 , but it can also be seen from Fig. [l] observing 
there only the trajectories starting inside [0, 1] are plotted). Since p is smooth, this 
implies that for R,T > there exists a positive constant cr^t such that 

p(t, x) > c R>T for (t, x) e [0, T] x [-R, R]. 



In particular, since ip e is given by (5.8) with r e small in C°, see (5.10), it follows 
that 

Cr,t 



(5.13) 



p s (t,x) > 



for (t,x) e [0,T] x [-R,R] 



for all e > sufficiently small (the smallness depending on T and R). 
that 

p 

X £ = u e (t, X s (t, x)), u £ = — , 



Recalling 



and that J e and p e are both smooth, it follows from (5.13) that u e is smooth as 
well inside [0,T] x [— R, R]. In particular, by the Cauchy-Lipschitz theorem, the 
Bohmian trajectories X s can never cross inside [0, T] x [-R, R]. Since by symmetry 
X £ (t, 1/2) = 1/2 for all t > 0, this implies in particular that, for all t G [0, T}\ 



X E (t,x)^l/2 Vx^l/2, 



X £ (t 1 x)^l/2 Vi^l/2. 

[0,T]x 



Letting e — >• we deduce that X £ X (locally) in measure on Q 
since otherwise the above property would give 

X(t,x)^l/2 Vx^l/2, X(t,x)^l/2 Vxs^l/2 

for allt ^ 0, which is not the case (see Fig. |IJ). This proves Assertion (iii). 



-R, R], 



□ 



Remark 5.5. 



Note that for \t\ > T* , i.e., after caustic onset, the Wigner measure 
t) \ ^froo- In particular, this shows that w is 



is given by ( |48| > for all (t, x) € (E t x 
insensitive to the Maslov phase shifts, since \a,j | 2 = \bj | 2 for all j = 1, . . . , N(t, x). 
The limiting Bohmian measure /3, however, incorporates these phase shifts in view 



20 



A. FIGALLI, C. KLEIN, P. MARKOWICH, AND C. SPARBER 



of the formula given in Theorem 4.2 However, as we have seen in Section 4.2 these 



phase shift do not enter in the classical limit of p £ and J £ . 

5.3. Extension to the non-zero potential case. In the case where V{x) ^ 
the situation becomes considerably more complicated, due to a lack of an explicit 



integral representation for the exact solution ip £ of ( 1.1 ). The only exception there 



from is the case V(x) = ±||x| 2 where one has Mehler's formula replacing (5.6) 



see, e.g., |llj . In order to proceed further in situations where V is a more general 
(sub-quadratic) potential, one needs to approximate the full Schrodinger propagator 

U e (t) = e~ iHt , with H £ = — A + V(x), 

for < e <C 1 by a semi-classical Fourier integral operator [T5] . Early results on this 
can be found in [131 118] , where the occurrence of caustics makes the approximation 
valid only locally in-time. This problem can be overcome, by considering a class of 
Fourier integral operators whose Schwartz kernel furnishes an e-oscillatory integral 
with complex phase and quadratic imaginary part, see |29l Theorem 2.1] for a precise 
definition. Using this, the authors of [25] construct a global in-time approximation 
of U £ (t) for potentials satisfying V E C£°(IR d ), i.e., smooth and bounded together 
with all derivatives (see also [211 ES] for closely related results with slightly different 
assumptions). By applying the stationary phase lemma to this type of (global) 
Fourier integral operator, one infers the following result, as a slight generalization 
of Theorem 5.1]: 

Fix a point (to,Xo) E (K* x K^) \ If, i.e., away from caustics, and as before 
denote by Yj(t, x) and j = 1, . . . , N = N(t, x) E N, the solutions of the equation 



x = X(t,y), where t t-t X(t,y) is the classical flow map induced by (1.6). Let 
{y 6 R d : \a (y)\ > 0}, be a sufficiently small neighborhood of 

(5.14) {n(*o,zo), • ■ • ,Y N (t Q ,x )} C R d } 

i.e., the points obtained by tracing back the classical trajectories intersecting in 
(t ,x ) £ (K t x R d ) \ ff. Then the solution of ( fTTT] > at t = t admits the following 
approximative behavior: 

N(t,x) 

(5.15) xjf{t ,x) e ^= + J2 a j (t ,x)e i * {m i^ x) - m ^ {to > x))/i e iS ^ to ' x ^ e +r e (t ,x), 

where the amplitudes Oj and the (real- valued) phases Sj are, respectively, given by 
(1.8) and (1.7) with Y replaced by Yj(t ,x), and mt(to,i) (resp. mj(t ,x)) is the 
number of positive (resp. negative) eigenvalues of the matrix V y X t (Yj(to, x)). In 
addition, the remainder r £ satisfies 

||r e (to,-)IU»(A) =0(e), 
where x E A C M. d is a sufficiently small neighborhood of xq e R d . The above 



result (the proof of which can be found in [5]) replaces Lemma 5.1 valid in the 
free case. Note however, that one only infers a local result in some sufficiently 
small neighborhood of xq € M. d , provided the initial amplitude oq is sufficiently 



concentrated on (5.14). In order to obtain an estimate for eVr e , we note that by 
applying the Hamiltonian H £ to (1.1), and having in mind that V E L°°(M. d ), we 



infer 

sup \\£ 2 &4> £ {t,-)\\ L 2 s; c, v< e 



where C > is independent of e. In view of (5.15), we consequently obtain that 



\ L 2 is uniformly bounded w.r.t. e and hence we can interpolate 

||eVr e ||£ 2 «SC||r £ || L2 ||e 2 Ar £ || i2 = 0(e), 
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to obtain ||eVr , e ||i2 = 0(y/e) — o(l), as required in Theorem 4.2 In order to apply 
the latter we also require VSj ^ VSfc for j ^ k 6 {1, . . . , N}. This follows, from 
similar arguments as has been done in the free case. Indeed, if the gradients were 
the same, by following backward the Hamiltonian flow we would get that the curves 
were starting from the same point, which is a contradiction. 



Thus, after using appropriate localization arguments, the multi-phase form ( 5.15 ) 
combined with Theorem |4.2| allows to infer the same qualitative picture for the 
classical limit of Bohmian dynamics in the case V ^ 0, as we showed above for the 
free case. Using the same notation as above, we can summarize our discussion as 
follows. 



Proposition 5.6. Let V € C£°(R ) and So satisfy Assumption 5.1 Let (t ,Xo) £ 
(R t x M~.) \ ^ , and assume that {y G M. d : \do(y)\ > 0} is a sufficiently small 
neighborhood o/{Yi(fo, xo), . . . , Yjsr(to, Xq)}- Then there exists a small neighborhood 
W C if x of (to, xq) such that j3 ^ w inside UxRp. In particular, the Bohmian 
trajectories (X £ , P £ ) do not converge locally in measure to the classical Hamiltonian 
flow. 

6. Numerical simulation of Bohmian trajectories 

In this section we shall numerically study the behavior of Bohmian trajectories, 
mainly in the regime < e « 1 and in particular in situations where caustics 
appear. Let us remark that the numerical implementation of Bohmian trajectories 
is used in applications of quantum chemistry, cf. [HI EI] • 

6.1. Description of the numerical method. For the numerical tracking of 



Bohmian trajectories (X 6 , P e ) it is necessary to solve the system (1.131 for a given 



solution tp E (t,x) of the Schrodinger equation (1.1 1. To this end, we will always 
consider initial data ipo G <S(R d ), i.e., rapidly decreasing functions. This allows to 
numerically approximate the solution tp E through a truncated Fourier series in the 
spatial coordinates by choosing the computational domain fi com sufficiently large, 
i.e., such that \tjj E \ is smaller than machine precision at the <5T2 com (we use double 
precision which is roughly equivalent to 10 -16 ). Thus the function can be period- 
ically continued as a smooth function with maximal numerical precision. In our 
numerical examples, we shall concentrate on the case of d = 1 spatial dimension. 
The x-dependence of ip e is consequently treated with a discrete Fourier transfor- 
mation realized via a Fast Fourier Transform (FFT) in Matlab. We thereby always 
choose the resolution large enough so that the modulus of the Fourier coefficients 
decreases to machine precision which is achieved in the studied examples for 2 10 to 
2 14 Fourier modes. This resolution enables high precision interpolation from x to 
X £ (see below). 

For the time-integration of the Schrodinger equation we shall rely on a time- 
splitting method. The basic idea underlying these splitting methods is the Trotter- 
Kato formula [3Hj, i.e., 

(6.1) lim ( e -tA/n e -tB/n\ n = e -t(A+B) 

where A and B are certain unbounded linear operators, for details see |27j . In 
particular this includes the cases studied by Bagrinovskii and Godunov in [3] and 



by Strang [37] ■ The formula (6.1) allows to solve an evolutionary equation 

d t u = (A + B)u, u\ t =o=u , 

in the following form 

U {t) = e ciAtA edl AtB e c 2 AtA e d 2 AtB . . . £ c k AtA^AtB^ 
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where (ci, . . . , c&) and (di, . . . , dk) are sets of real numbers that represent fractional 



time steps. In the numerical treatment of (1.1) we shall use a second order Strang 
splitting, i.e., Cj = di = 1 for all i except for c\ = dk = 1/2- The Schrodinger 
equation is consequently split into the following system: 

£ 2 

iedtu + —d xx u = 0, iedtu — V(x)u. 

The first equation can then be explicitly integrated in Fourier space, using two 
FFT's. The second equation can explicitly be solved (in physical space) in the 
form 

u(t,x) =e- ttvi - x)/s u . 



Next, in order to solve the Bohmian equations of motion (1.13) for a given 
tp £ (t,x), we need to interpolate between the coordinate x, in which ip e is given, 
and the coordinate X s . For this we use that the x-dependence of V> £ is treated by 
Fourier spectral methods. Thus we can apply the representation of ip £ in terms of 
truncated Fourier series not only at the collocation points for which the formulae 
for the discrete Fourier transform hold, but at general intermediate points. The 
main drawback is that for such points there is no FFT algorithm known and the 
transformation is thus computationally more expensive. But since we only need to 
track a limited number of trajectories X s and since this interpolation method is 
of high accuracy, our approach is more efficient than, say, a low order polynomial 
{spline) interpolation (as used, e.g., in [14]). In order to obtain the Bohmian 
momentum P £ we interpolate, x <-> X s within tp e (t,x) and d x ip e (t,x), for fixed 
time t G M. To this end, we note that the latter is of course determined in Fourier 
space. We consequently compute P e through 

We test the accuracy of the interpolation by comparing different numbers of Fourier 
modes for the solution of the Schrodinger equation for a given set of computed 
trajectories. Once machine precision is assured for ip e (i.e., the modulus of the 
Fourier coefficients decreases below 10~ 12 , in our case), the difference between 
different interpolates can be shown to be of the same order. Thus we can conclude 
that the spatial resolution of the trajectories is of the order of 10~ 12 , much better 
than plotting accuracy. 



The time integration of the first equation of the system (1.13) is performed 
with an explicit scheme (here, we shall use a standard fourth order Runge-Kutta 
method). This allows to compute the right-hand side of this equation with the 
already known values for X s at the previous time step. Note that we compute the 
solution to the Schrodinger equation either exactly in time (if V(x) = 0) or with 
second order time splitting for each stage of the Runge-Kutta scheme (whenever 
V{x) 7^ 0). We shall test the accuracy of the time integration scheme by assuring 
that the difference of the numerical solution for N t time steps to the solution for 
2N t time steps is smaller that 10 -4 and thus much smaller than plotting accuracy. 
Typically we use Nt — 10 4 . In addition the accuracy of the splitting scheme is 
tested as in [28] by tracing the numerically computed energy £^ um (i) which due 
to unavoidable numerical errors is indeed a function of time. In our examples, the 
relative conservation of £ l „ um (t) is ensured to better than 10 -7 implying again an 
accuracy of more than 10~ 5 . 

Remark 6.1. For efficiency reasons, the computation of the trajectories X 6 is 
done at the same time for all X s . Thus, in principle, it could happen that the 
identification of the trajectories in the examples below do not reflect the actual 
dynamics. By tracing also individual trajectories, i.e., by computing just one X s 
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per run, we nevertheless are able to ensure that this is not the case and that the 
shown trajectories are indeed the correct ones. In particular our numerical code 
captures the physically imporant property that Bohmian trajectories do not cross, 
see, e.g., [14] (see also the proof of Theorem 5.4 (hi)). This is indeed a delicate 



issue in other numerical approaches where the system ( 1.6 ) is numerically integrated 
with ( 1.3 1 and ( |1.4[ ) instead of ( |1.1[ ), and where different interpolation techniques are 
used. The latter have to be chosen in a way to avoid the crossing of the trajectories 



(see Section 6.2.1 below) 



6.2. Case studies. In the following we shall illustrate our analytical results by 
numerical examples, starting with the (globally smooth) case of semiclassical wave 
packets, which has already been treated in an earlier paper [33]. We shall then also 
consider the case of i/jq exhibiting vortices before we finally deal with WKB initial 
data producing caustics. 

6.2.1. Vortices. Before studying the regime < £ < 1 we shall show that our 
numerics displays an important non-crossing property of Bohmian trajectories. In- 



deed, it is well known that solutions to the Schrodinger equation ( 1.1 ) in general can 
have nodes, i.e., points at which the wave function vanishes. Due to the superfluid 
property of i> e such nodes represent quantum mechanical vortices. At such points, 
the Bohmian trajectories X e are not well defined, but since P £ does vanish as well 
at these points, there is a natural analytic continuation of the trajectories through 
such nodes. In the following, we shall numerically study the example given in [6]. 
More precisely, tp e is given by the superposition of the ground state and the second 
excited state of the harmonic oscillator (we also put e = 1 in this example), i.e., 

%l>{t,x) = (1 + (1- 2x 2 )e- 2lt ) e~ x2 - u/2 . 

This wave function vanishes for x = and for all times t — (2k + l)n/2, with k E Z. 
To treat the limit '0/0' numerically, we add some quantity of the order of the 
rounding error to the wave function which will consequently provide the limit with 
an error of the order of the unavoidable numerical error. The resulting trajectories 
can be seen in Fig. [2] Note that indeed, all trajectories avoid the vortices at t = n/2 
and t — 3tt/2, only the trajectory for x — passes through these nodes. 




Figure 2. Bohmian trajectories for e = 1 in a harmonic oscillator 
potential V(x) = \x 2 with tp given as a superposition of the ground 
state and the second excited state. 
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6.2.2. Semiclassical wave packets. In |33) a result similar to Theorem 3.1 (ii) is 
proved, for the case of semiclassical wave-packets (see also |16) for a closely related 
study). The corresponding initial data are of the form 

(6.2) %(x) = £- d / 4 a f^p) e ik < x - x °V e , a e S(R d ; C). 

The main differences between WKB states and semiclassical wave packets are that 
for the latter, the particle density concentrates in a point, i.e., 

Pl{x) e ^ S(x-x ), inV'(R d ), 

and that the corresponding semiclassical approximation does not exhibit caustics, 
cf. |12j for more details. This in particular implies that for semiclassical wave 
packets one can prove convergence of the Bohmian trajectories on any finite time- 
interval [33]. An example for such a situation (with fco = 0) can be seen in Fig. [3j 
The corresponding classical trajectories would be just lines parallel to the i-axis. 



Figure 3. Bohmian trajectories for wave packet initial data of the 
form (|eT2j> with k Q = 0, x = 1/2, a Q (z) = eT** and e = 1CT 3 . 

Since these data do not lead to a caustic, there is just a slight defocusing effect to 
be seen with respect to the classical trajectories. 

6.2.3. Caustics. In this last subsection we shall, finally, present examples exhibiting 
caustics in the classical limit. To this end, we shall first study the case where the 
caustic is just one single point, i.e., a situation in which all classical trajectories 
X(t,y) cross at (x*,T*) el ( x M. x . As a particular example, we shall consider the 
harmonic oscillator with potential 

(6.3) V{x)= i ^ x _ i y 

and an initial data in the form 

(6.4) ^(x) = 6-25(^-1/2)^ 

i.e., a WKB state with Gaussian amplitude and Sq(x) = 0. Then, the classical 
trajectories X(t,y) all intersect in one point as can be seen in Fig.|4j 

The same situation for the Bohmian trajectories X £ (t,y) and e = 10~ 3 can be 
seen in Fig. [5j The closeup of the region of intersection when e = clearly shows 
that the trajectories come close to x* , but keep a finite distance from it except for 
the one trajectory which is parallel to the t-axis and goes straight through x* . The 
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Figure 5. Left: Bohmian trajectories X e {t,y) for the harmonic 
oscillator potential (6.3) and V>o given by ( |6.4[ ). Right: A closeup 
of the central region near x* . 



solution "0 6 is periodic in time and shows a breather- type behavior with a large \ip 6 \ 
at the caustic. We show only a half-period of this periodic motion. 
Next, we consider the case V(x) = with WKB initial data 

(6.5) i/fo(x) = e -25(z-i/2) 2 e *So(z)/ £j 5 (x) = -ilncoshf5a;- 

as in [31] . i.e., the same amplitude as before but with nonzero initial phase. The 
time dependence of the density p shows a strong maximum followed by a zone of 
oscillation inside a break-up zone as can be seen in Fig. [6} In this case, the classical 
trajectories X(t, y) will lead to a diffuse caustic as depicted in Fig. [I] For finite 
e, the Bohmian trajectories X e (t,y) obviously do not cross, but there are rapid 
oscillations within the caustic region as can be seen in Fig. [7J 

However, oscillations do not only appear in the trajectories, but also in the 
momentum P £ (t,y) = u £ (t, X s (t,y)) along any trajectory X e which is "deflected" 
at the caustic region. This can be clearly seen in Fig. [8] where several P £ are plotted 
along the corresponding trajectories X £ . The oscillations within P £ are reminiscent 
of so-called dispersive shocks, as observed, e.g., in the Korteweg-de Vries equation 
with small dispersion, see for instance (52] and references therein. This is even 
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Figure 6. Density p = \ip e \ 2 for V(x) = and ip^ given by (6.5) 





Figure 7. Left: Bohmian trajectories X e (t,y) for V{x) = and 
ipQ given by (6.5 1. Right: A closeup of the central region. 



more visible in Fig. [9] where the oscillations on the left most trajectory in Fig. [8] 
are shown in dependence of t, thus in a projection onto the i-axis. 
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